Large-scale investigations of Neolithic settlement dynamics in Central Germany based on machine learning analysis: A case study from the Weiße Elster river catchment

The paper investigates potentials and challenges during the interpretation of prehistoric settlement dynamics based on large archaeological datasets. Exemplarily, this is carried out using a database of 1365 Neolithic sites in the Weiße Elster river catchment in Central Germany located between the southernmost part of the Northern German Plain and the Central Uplands. The recorded sites are systematically pre-processed with regard to their chronology, functional interpretation and spatial delineation. The quality of the dataset is reviewed by analyzing site distributions with respect to field surveys and modern land use. The Random Forests machine learning algorithm is used to examine the impact of terrain covariates on the depth of sites and pottery preservation. Neolithic settlement dynamics are studied using Site Exploitation Territories, and site frequencies per century are used to compare the intensity of land use with adjacent landscapes. The results show that the main trends of the Neolithic settlement dynamics can be derived from the dataset. However, Random Forests analyses indicate poor pottery preservation in the Central Uplands and a superimposition of Neolithic sites in the southernmost part of the Northern German Plain. Throughout the Neolithic the margins between soils on loess and the Weiße Elster floodplain were continuously settled, whereas only Early and Late Neolithic land use also extended into the Central Uplands. These settlement patterns are reflected in the results of the Site Exploitation Territories analyses and explained with environmental economic factors. Similar with adjacent landscapes the Middle Neolithic site frequency is lower compared to earlier and later periods.


Introduction
Humans always depended on natural resources and environmental conditions and altered Holocene landscapes themselves by different activities such as settlement and infrastructure construction, fire ignition, tilling, herding or mining on various time scales [1][2][3][4]. Therefore, the diachronic reconstruction of former human activity in a region is crucial to understand former human-environmental interactions. There are two different methodological approaches to reconstruct former human impact in a region: (i) Using geoscientific methods, regional human impact is often derived from sediment archives using paleoenvironmental proxies such as pollen, biomarker, charcoal abundance or gastropods [5][6][7][8]. However, although such studies often show a high temporal resolution, upscaling is difficult and they hardly give information about spatial patterns of former human activity. (ii) Spatially resolved information can be obtained by the systematic study of regional settlement patterns based on archaeological data [9][10][11][12]. In Germany, a systematic recording of prehistoric sites in inventories, catalogues and maps started at the beginning of the 20th century against the background of an increasing institutionalisation of archaeological research and heritage management [13][14][15]. This enabled archaeologists to conduct large-scale diachronic studies, i.e. to investigate the prehistoric settlement dynamics in various landscapes by comparing maps with site distributions from different time periods [16,17]. Most of these studies used archaeological and geographic data to discuss long-term changes in the interactions between early farming societies and their environments. Usually, this was achieved by plotting archaeological sites against different types of soils and land cover [18][19][20][21][22][23][24][25][26][27][28].
In principle, this research approach is still being pursued to the present day. However, the methods applied have been improved: (1) additional terrain covariates are used to investigate settlement patterns. Until today, the most frequently studied terrain covariates include elevation, relief position, slope, exposition, distance to rivers, precipitation, temperature and soil varieties [12,29]. Occasionally, soil fertility, usable field capacity as well as cation exchange capacity and the suitability of a soil for ploughing are taken into account as well [30][31][32]. (2) Furthermore, there was a shift towards statistical analyses to improve the reproducibility of research results [33][34][35][36][37]. (3) In addition, standard methods have been established to enable an evaluation of the representativity of archaeological maps. These include a semi-quantitative review of the local research history, a discussion of site distributions with respect to modern land use, the nature and intensity of archaeological field work as well as statistical investigations on the impact of erosion, colluviation and weathering conditions on the preservation of prehistoric sites [38][39][40]. Furthermore, comparisons of archaeological site distributions to random point distributions as well as (predictive) modelling approaches are applied to the study of former settlement dynamics [11,37,[41][42][43]. (4) Finally, there was a shift from spatial analyses based on point coordinates to studies based either on Site-Catchment Analysis or Site Exploitation Territories (SET) [12,[44][45][46][47][48][49][50][51].
Up to now, a large number of diachronic studies dealing with long-term human-environment interactions between the Neolithic and the Middle Ages has been carried out (Fig 1). However, these are distributed unevenly across Germany, with a focus on landscapes along the rivers Danube [10,12,52,53] and Main [9,29,54, 55] as well as in Central Germany [33,41,56,57]. So far, only one such study has been carried out in the Northern German Plain [58]. Aside from their uneven spatial distribution, it is striking that most studies hardly refer to each other, e. g. regional results are rarely discussed in the context of supra-regional developments [10,12,29,32].
The present study was carried out in the frame of an interdisciplinary geoarchaeological research project that systematically investigates the interplay between the Holocene geomorphological floodplain and slope dynamics, variations of the prehistoric and medieval settlement dynamics and climate changes in large parts of the loess-covered catchment of the Weiße Elster river in Central Germany [64,65]. We present first results from the large-scale prehistoric to medieval database that was built up during this project, and discuss potentials and challenges • assessing the current state of knowledge on the Neolithic chronology and culture sequence in the study area by reviewing the local history of (field) research • evaluating Early, Middle and Late Neolithic site distributions with regard to archaeological field surveys and modern land use • analysing the interrelation between several terrain covariates, and the depth of Neolithic sites as well as local conditions for the preservation of Neolithic pottery and stone tools, to assess possible influences of geomorphological terrain dynamics and weathering on the archaeological record • identifying possible biases in the archaeological record • analysing Early, Middle and Late Neolithic settlement dynamics by using site densities, SETs and site frequencies • discussing the results in the light of observations from earlier studies in Central Germany

Study area
The study area is situated in Central Germany along the borders of the federal states Saxony, Saxony-Anhalt and Thuringia (Fig 2). It is defined by the immediate catchment of the Weiße Elster river between the city of Leipzig in the north, and the German-Czech border in the Elster Mountains (German: Elstergebirge) in the south. It covers an area of approx. 3026 km 2 . Geographically, it is divided into two distinctive landscapes: the northern part belonging to the southernmost part of the Northern German Plain (German: Norddeutsches Tiefland), and the southern part in the Central Uplands (German: Mittelgebirgsschwelle) (S1 Fig) [66][67][68][69]. Geologically, the southern part of the Northern German Plain correlates with Cenozoic deposits, whereas Mesozoic to Palaeozoic formations form the southern Central Uplands [70][71][72]. In the northern part, the average elevation ranges between 100 and 200 m above sea level (a.s.l.), and gently rolling slopes, fertile chernozem/phaeozem and luvisol soils on loess are found [73][74][75]. The mean annual air temperatures range around 10˚C and the mean annual precipitation range between 550 and 650 mm [76,77]. In the southern Central Uplands, the elevation rises up to 700 m a.s.l., resulting in lower mean annual air temperatures (9-7˚C), increasing amounts of annual precipitation (up to 1050 mm) as well as longer periods of winter and frost [76][77][78][79]. Here, the topography is more accentuated with deeply cut river valleys and steep slopes. Furthermore, the landscape upstream is less favourable for agricultural use due to lowyielding dystric cambisols and stagnic gleysols towards higher altitudes [80,81].
To analyse the spatial distribution of the recorded sites, we distinguished between sites with an accurate and a symbolic location. The first category applies to sites that can be located in the field. Symbolic coordinates were assigned to sites that cannot exactly be located due to missing or poor geographic information. In these cases, either the presumed location or the church of the village in whose district the site was discovered was taken as the point of reference [99].
Based on the information given in the local area files and literature, three types of sites were differentiated: (I) Settlements are indicated by pits, post-holes, grinding stones, loom weights, spindle whorls or pottery fragments. (II) Burial sites were defined by the presence of human Study area with respect to borders of federal states, international (intl.) borders, local topography and large urban areas mentioned in the text. The outlines of the study area were extracted from the Catchment Characterisation Model database [82]. The political borders are drawn according to the Saxon State Office for Environment, Agriculture and Geology, the Saxony-Anhalt State Office for Surveying and Geoinformation and the Thuringian State Office for Land Management and Geoinformation. Rivers and urban areas are derived from the urban morphological zones dataset and the European Catchments and Rivers network system [83] remains. (III) Stray finds included single objects or assemblages with little or no context information. Accordingly, these finds may be either lost objects or material remains from settlements or burial sites that have not yet been identified as such [58,100]. Therefore, a site was registered twice in our database whenever a burial was documented in a settlement.
Due to frequent construction measures and/or archaeological field work in a relatively small area, the residuals of a settlement or burial site may be discovered in multiple spots close to each other. If each of these archaeological observations were taken as hints to an independent site, the erroneous impression of a densely settled area might emerge, when in fact only the intensity of archaeological observations is reflected. Therefore, we aggregated sites that (I) were less than 250 meters apart from each other, (II) associated with the same or a compatible dating (e.g. "Neolithic" + "Neolithic", Early Neolithic + "Neolithic", etc.), and (III) interpreted in the same or a compatible manner (e.g. settlement + settlement, settlement + stray finds, etc.). Since the actual spatial extent of prehistoric sites remains unknown until they have been excavated individually, it is important to acknowledge that this aggregation is ultimately based

PLOS ONE
Large-scale investigations of Neolithic settlement dynamics based on machine learning analysis on an arbitrarily chosen standard. Nonetheless, without this measure analyses of site distribution patterns are prone to result in an overrating of geographical parameters [10,32,101,102].
3.1.3 Pre-processing of spatial data. With respect to the spatial resolution of the archaeological data, continuous raster data and categorical grids derived from vector data were resampled to a resolution of 250×250 m (S1 Table). In addition, corrections of the applied digital elevation model (DEM) [81] were necessary, because the topography south of Leipzig has been altered by modern open-cast mines [103,104]. Therefore, we georeferenced pre-mining topographic maps (S2 Table) and digitized their contour lines. Thus, we were able to remove remnants from open-cast-mines such as pits, fills or mounds in an area of approx. 640 km 2 (S2 Fig). In addition, the data for water bodies and rivers were edited based on the georeferenced older topographic maps. This was necessary since several river courses were altered as a result of the open-cast-mines, and because the datasets also contained artificial water bodies that cannot be used for an analysis of the prehistoric settlement dynamics (S3 Fig). The result is an approximation to prehistoric and medieval conditions, as the exact course of each river is subject to continuous changes [105,106].

Archaeological source criticism.
To evaluate Neolithic site distribution patterns and thus the reliability of the results that can be derived from them, an analysis of the archaeological data is necessary. In principle, it must be assumed that the distribution of prehistoric sites is influenced by (I) archaeological filters, i.e. the type and intensity of archaeological field work, (II) geographical filters such as modern land use, erosion, colluviation or weathering and (III) culturally intrinsic filters, i.e. the actual understanding of land use strategies in prehistory. The aim of archaeological source criticism is therefore to approximate the degree of distortion caused by the first two filters, and thus to evaluate the overall quality of the archaeological database [86,[107][108][109].

State of local research.
The local state of archaeological research will be presented in the form of key data describing the basic composition of the Neolithic record according to our database. These include the number of recorded sites with respect to their dating, the interpretation of their individual function, and the accuracy of their location as well as their status of publication [9][10][11][12]32, 53-59].
3.2.2 Intentionality of site discoveries. It is possible to distinguish between intentional and non-intentional modes of site discoveries: The first category includes archaeologically motivated activities such as excavations, field surveys or aerial photography. In the second category, activities without archaeological motivation are subsumed, including construction measures, the extraction of raw materials (mining) or activities related to agriculture and forestry [110]. The discrimination between these categories is important, because non-intentional modes of discovery represent a random sample with regard to the spatial and chronological distribution of prehistoric sites. Accordingly, the analysis of the intentionality is crucial for the understanding of the representativeness of the dataset [9].
Furthermore, territories of individual collectors might affect site distribution patterns. Therefore, it is necessary to discuss the Neolithic site distribution with respect to areas that were investigated by field surveys or aerial photography [10,12,32,39]. Taking into account individuals who discovered more than nine sites [100], we modelled the territories covered by each method separately from each other on the basis of the largest empty circle [111][112][113]. This discrimination was done because aerial photography may result in the discovery of archaeological features on or slightly below the recent surface. However, in contrast to field surveys, this method cannot contribute to the discovery of small stone tools and pottery fragments. For each approach, the investigated territories were described by 1.5 and 2.5 km isolines. It should be borne in mind that the modelled territories were derived from data of surveys that successfully resulted in site discoveries. Accordingly, the modelled territories are not necessarily congruent with the areas actually surveyed, but rather smaller.

Geographical filters
3.3.1 Site distribution in relation to archaeological site visibility and preservation of material remains. Local geomorphic terrain dynamics are known to affect archaeological site visibility and the preservation of material remains: Firstly, colluviation and alluviation can lead to the superimposition and thus the protection of prehistoric sites at the cost of their visibility [38, 114,115]. This can be discussed by using information about the level of an archaeological site below the recent surface that can be drawn from the activity leading to its discovery [9,39,58,59]. Secondly, erosion may lead to the uncovering of prehistoric sites, resulting in an improved accessibility at the cost of reduced preservation conditions or even their complete destruction due to an increased exposure to weathering and ploughing [107,108]. To examine this aspect, we differentiated between sites where pottery (and stone tools) were recorded, and those where only stone tools were discovered.
To discuss possible effects of terrain covariates (S1 Table) on the depth of Neolithic sites as well as the presence/absence of pottery, we used the Random Forests (RF) modelling algorithm, i.e. an ensemble learning method for classification and regression [116]. RF is suited for spatial archaeological research and predictive modelling, because it can handle complex spatial and contextual datasets with weak, imbalanced or noisy inputs. RF consists of a large number of individual decision trees that form an ensemble, where each individual tree is based on a distinct set of binary decision rules. In addition, RF reduces over-fitting by averaging (regression) or aggregating by majority vote (classification) of all individual tree results [117,118]. A description of the RF algorithm is provided below: i. Initially, a random subset of the original dataset was generated by using the bootstrap method (random sampling with replacement) for each individual tree. Thus, we did not use the entire dataset for learning, but the same number of samples. Usually, the original dataset is represented in the sample set with~63% of unique and~37% of redundant samples. Thus, the remaining samples not represented in the training set are called out-of-bag (OOB) samples and used for validation. ii. A second-level random component was used to build the individual tree. Therefore, at each node the best binary split on the training data was performed, based on a random subset of predictors (S1 Table). The number of covariates (mtree) can be chosen by the user. However, Breiman [116] suggests to use the square root of the number of predictors. Based on a binary recursive partitioning method a distinct set of binary rules is generated within all individual trees, with the aim to separate each single class (classification case) in a terminal node. This is done until no further partitioning is possible, e.g. a minimum set of samples is represented in a node or the node contains only samples with one single class label (perfect split).
iii. If used for a classification, the results from each decision tree were aggregated on the basis of a majority vote.
The algorithm provides a distinct set of accuracy measurements to evaluate the performance of the RF model [119]. These are based on the confusion matrix (Fig 4) [120], which consists of four specific conditions (true positive, false positive, false negative, true negative) related to the predicted OOB output. This allows a detailed analysis of the classification accuracy. The out-of-bag error rate is calculated as follows: OOB error rate ¼ sum of misÀ classified observations sum of observations in training set ð1Þ The Kappa value measures the agreement between predictions and observed cases. This is achieved by comparing the overall accuracy to an expected random chance accuracy. Therefore, the Kappa score is particularly useful on classifications with imbalanced data [121][122][123][124]. The results range from -1 to 1 with the latter being a perfect classification, followed by strong agreement (>0.8), substantial agreement (between 0.4 and 0.8) and poor agreement (<0.4). A score of 0 indicates that the classification is as good as random. Values below 0 point to no effective agreement or even disagreement [125,126].
Cohen's Kappa is calculated using the following formula [124]: Where Pr(a) represents the actual observed agreement (true positives) and Pr(e) chance agreement.
Pr e ð Þ ¼ Where cm 1 represents column 1 marginal, cm 2 column 2 marginal, rm 1 row 1 marginal, rm 2 row 2 marginal and n the number of observations. In addition to the described measures of agreement in a categorical use-case (classification), we used additional validation measures by deriving so-called measures of effectiveness such as recall, precision and the F 1 -measure based on the confusion matrix, which are important measures for information retrieval [120,127].
Recall (Rc) describes the relation between positive (true positive: TP) and negative (false negative: FN) predicted class values, and thus the probability that a pre-classified sample is actually predicted, i.e. it represents underestimation. In a prediction without erroneous assignments (FN), the value is 1. The recall is calculated as follows: In contrast, precision (Pc) describes the probability that an estimated prediction result is actually also pre-classified within the training dataset (true negative: TN), i.e. it represents overestimation. In a prediction without erroneous assignments (FP), the value is 1. The precision is calculated as follows: To quantify the overall model quality in a composite way the F 1 -measure is frequently used. The F 1 -measure represents the harmonic mean between under-and overestimation, and is calculated as follows: Additionally, RF provides a measure of importance of each covariate by calculating their mean decrease accuracy [Formula 6], i.e. the impact of each variable on the prediction error when it is not taken into account [127,128]. This enables in-depth analysis of interaction from the covariate space related to our classification problem. Breiman [116] suggests the difference in prediction accuracy before and after permuting a single covariate, averaged over all trees, as a measure for variable importance. The mean decrease of accuracy is calculated as follows [128]: where � B ðtÞ is the out-of-bag sample for a tree, y ðtÞ i ¼ f ðtÞ ðx i Þ the predicted class for observation i before, and y ðtÞ i;p j ¼ f ðtÞ ðx i;p j Þ the predicted class for observation i after permuting its value of variable X j, i.e. with x i;p j ¼ ðx i;1 ; . . . ; x i;jÀ 1; x p j ðiÞ;j ; x i;jþ1 ; . . . ; x i;p Þ: Following Breiman [116], we used the square root of the number of predictors (mtree) and a maximum number of individual trees of 1500. Even though the general precondition of a RF algorithm is robust against multi-collinearity [116,129,130], we only used covariates with correlation of less than 0.8 (S1 Table; S4 Fig). To avoid statistical coincidence, all analyses were repeated ten times to calculate average values for the OOB error rate, Kappa, Precision, Recall and F 1 -measures. The statistical analyses were carried out using the packages randomForest [131] and caret [132,133] in R version 3.6.3 [134].

Site distribution in relation to modern land use.
Various large-scale studies showed that modern land use has an impact on the accessibility and preservation of prehistoric sites [10-12, 32, 39, 42, 135-137]. In recent years, CORINE Land Cover (CLC) raster data have been successfully used to evaluate the relationship between archaeological site distribution patterns and modern land use [10][11][12]42]. Based on satellite imagery, the dataset originally discerns 44 classes of modern land use with a spatial r esolution of 100×100 meters [83]. However, it is recommended to aggregate these classes for archaeological purposes [10,12,100]. Here, we distinguished seven types of land use: urban areas, forests, arable land, grassland, water bodies, bogs or swamps, and landfills or dumpsites.
We applied the Chi-squared test to assess the interrelation between the relative frequency of land use classes and Neolithic sites. This was achieved by calculating the number of sites expected to be affiliated with each land use class, and then comparing these values with the observed distribution of Neolithic sites over these classes. The significance of the deviations between the expected and observed sites was described by the χ 2 value, calculated as follows: Where k is the number of categories, O i the observed number of cases in category i, and E i the expected number of cases in category i. The sum of all χ 2 values was compared with the critical χ 2 value p [138][139][140].

Site frequency.
Due to the fact that the duration of each Neolithic period varies, we used the site frequency per century (German: Fundstellenfrequenz) to enable a comparison of the settlement intensity between the Neolithic periods. It was calculated as follows [34, 48]: Where S n is the number of sites and P n is the duration of the respective period in years. The site frequency can be used to compare varying settlement intensities on a supraregional scale [12,39,141]. However, the comparability between individual study areas is limited due to their varying size, which in turn affects the number of sites that can be recorded.
To improve the comparability of study areas of different sizes, we calculated the average site frequency for an area of 250 km 2 . The site frequency was processed as follows: Average site frequency per 250 km 2 ¼ Fr Where Fr is the local site frequency and Area the spatial extent of the respective study area in km 2 .

Site distribution and density.
To describe the spatial distribution and density of the recorded Neolithic sites, we calculated 1.5 and 2.5 km isolines based on the largest empty circle between sites with accurate coordinates [111][112][113]. This approach allowed the identification of core areas of Neolithic settlement on an empirical basis for each period. By comparing the extent and distribution of these areas, the settlement dynamics can be derived [142,143].

Site Exploitation Territories (SET) and settlement dynamics.
We applied the concept of SET to investigate the Neolithic settlement dynamics. This concept was developed to study archaeological sites in the context of their environment [144][145][146]. The term refers to a time-distance based territory that is presumably visited for daily subsistence [147][148][149]. Accordingly, the local topography (slope) influences the shape of SETs. In landscapes with gently rolling slopes SETs tend to circular shapes, whereas due irregular forms are more common in mountainous regions [113,150,151].
SETs are based on the surmise that each site has an optimal geographic location with respect to its economic function. Consequently, it is expected that mobile groups whose subsistence was pasture farming preferred locations favourable for grazing, whereas settlements from sedentary societies are expected to be located in areas suitable for agriculture [147,149,152]. With respect to anthropological field studies, it is assumed that the intensity of resource and land exploitation decreases with increasing distance from the settlement. In the case of sedentary societies, the territory used on a daily subsistence basis is assumed to be in the area that can be reached within one hour's walking distance [147][148][149]152]. However, the most important activities related to arable farming are considered to take place in the immediate vicinity of a settlement, i.e. in an area that can be reached within 10 to 15 minutes walking [148,149,[153][154][155][156]. In our study area this was confirmed by investigations of Neolithic wells from Eythra and Droßdorf. Due to good preservation of organic material, remains of cultivated plants, field weeds and insects were discovered in the wells that indicate agricultural activities in the direct vicinity of the settlements [157][158][159].
We modelled for all Neolithic sites with accurate coordinates their individual SET that can be reached within ten minutes [160]. The script implements the hiking function developed by Waldo Tobler [161] which has been used before in archaeological studies [162,163]. This function is an empirical model based on marching data of the Swiss military, taking into account multiple factors such as vegetation, individual physical fitness, length and quality of the path, altitudinal difference, weather conditions, darkness as well as marching competence and luggage [164].
Here, SET modelling was based on processed SRTM DEM data with a resolution of 250×250 meters, an estimated walking distance of 5 km/h on flat terrain as suggested by Tobler [161], and a damping of walking speed on slopes with an inclination of >15˚ [12]. The settlement dynamics were investigated by comparing two sets of grid statistics for each period: one derived from all sites dating to the respective period, and the other derived exclusively from settlements.
To further characterize the SETs, the following terrain covariates were used: area covered by the SET, elevation above sea level, slope and valley depth (i.e. height above nearest river) as derived from the processed DEM, river distance derived from the processed data and percentage of soils on loess (S1 Table). All these terrain covariates are subject to both natural and anthropogenic changes [105,106,[165][166][167]. In an ideal scenario, GIS analyses of prehistoric site distributions are based on individual datasets that are representative for each period. However, as there are no datasets that meet this standard, archaeological studies have to rely on available data [10,12,32]. With regard to the preceding remarks, the results of the SET analyses are to be understood as results of model calculations. They represent an approximation to archaeological data, and provide a basis for the discussion of settlement dynamics.

State of local research.
Following the systematic processing of all available information from the literature and the local area files, our record of the Neolithic in the study area consists of 1365 sites, including 361 Early, 82 Middle and 298 Late Neolithic sites (S3 Table). Due to the rather unspecific nature of the retrieved artifacts, the remaining 624 sites can only be described as "Neolithic". The Early Neolithic is represented by 208 settlements, 8 burial sites and 145 stray finds. In contrast, 47 settlements, 19 burial sites and 16 stray finds are associated with the Middle Neolithic. There are 74 settlements, 96 burial sites and 128 stray finds dating to the Late Neolithic. In total, 769 Neolithic sites are exactly located, including 340 settlements, 123 burial sites and 306 stray finds. Symbolic coordinates were assigned to 70 settlements, 26 burial sites and 500 stray finds (S3 Table). Not more than 703 sites have been published until today, which underlines the paramount value of using local area files in largescale studies (S4 Table).

Intentionality of site discoveries.
The mode of discovery is documented for 843 sites (S5 Table). Their majority is associated with non-intentional modes of discovery (n = 527). Construction measures (n = 201) and the extraction of raw materials (n = 182) take the largest share of this group. In contrast, 313 sites are linked to intentional modes of discovery, mainly field surveys (n = 288).
Based on the archaeological database, 12 individuals were identified who repeatedly carried out field surveys in the study area and discovered at least nine sites. In total, these individuals registered 263 prehistoric and medieval sites that were used to model the surveyed territories ( Fig 5). In total, 125 Neolithic sites were used for modelling the surveyed territories. This represents 30% of all Neolithic sites registered in those territories (n = 422), and 69% of all Neolithic sites discovered on the surface in these areas (n = 182). Accordingly, the field surveys resulted in an increase in the density of Neolithic sites, especially between Gera and Zeitz ( Fig 5).
Finally, territories for aerial surveys were calculated on the basis of 74 sites south of Leipzig and in the vicinity of Plauen. However, the aerial surveys did not affect Neolithic site distributions ( Fig 5).  Table). See Fig 2 for Table 1).
The majority of buried sites is located in the northern study area in the southernmost part of the Northern German Plain, while sites on the surface are more frequent in the Central Uplands (Fig 6A). A comparison of the material groups shows a concentration of sites with pottery in the northern study area. In contrast, in the area around Gera and further south in the Central Uplands there are almost exclusively sites without pottery (Fig 6B). The vast majority of discovered Neolithic pottery in the northern study area was buried at the time of its discovery. In contrast, in the southern Central Uplands Neolithic pottery tended to be on the surface (Fig 6C). The distribution of sites with stone tools but without pottery shows a concentration in the vicinity of Gera, and extends further south into the Central Uplands. A contrast between buried sites and discoveries on the surface is not recognisable (Fig 6D).  Table). See Fig 2 for a location Table). Therefore, based on the selected terrain covariates (S1 Table) the RF algorithm could reproduce the distribution pattern of the depth of sites (Fig 6A). According to the mean decrease accuracy, the five most important terrain covariates are (1) elevation, (2) copper distribution, (3) the proportion of coarse fragments and (4) sand in the soils as well as (5) the slope length and steepness factor (LS-factor).
In addition, the RF algorithm was able to reproduce the presence-absence-pattern of Neolithic pottery (Fig 6B), and resulted in an average OOB estimate of the error rate of 22.64%, Kappa score of 0.53, Precision of 0.79, Recall of 0.81 and F 1 of 0.79 (S7 Table). The five most important terrain covariates were (1) elevation, (2) copper distribution, (3) LS-factor, (4) slope and (5) the proportion of coarse fragments in soils.
Furthermore, the RF algorithm was able to reproduce the distributions of pottery with regard to site depth, and resulted in an average OOB estimate of the error rate of 19.40%, Kappa score of 0.55, Precision of 0.75, Recall of 0.63 and F 1 of 0.68 (S8 Table). Hence, based on the selected terrain covariates the RF algorithm reproduced the classification of buried and non-buried pottery (Fig 6C). The five most important terrain covariates were (1) copper distribution, (2) the proportion of coarse fragments in soils, (3) elevation, the (4) wind erodible fraction of the soil and (5) slope.
With regard to the distribution of (non-)buried stone tools, the RF model resulted in an average OOB estimate of the error rate of 24.54%, Kappa score of 0.21, Precision of 0.75, Recall of 0.89 and F 1 of 0.81 (S9 Table). Accordingly, the low Kappa score demonstrates that unlike the preceding analyses the RF algorithm could not reproduce the classification of (non-)buried stone tools (Fig 6D) based on the selected terrain covariates.

Site distribution in relation to modern land use.
The recorded Neolithic sites are unevenly distributed over modern land use classes (Fig 7). According to the Chi-squared test this observation is statistically significant ( Table 2). This is mainly due to the fact that approx. twice as many sites as expected were discovered in urban areas and mining areas. Also on arable land the number of recorded sites exceeds the number of expected sites. In contrast, some land use classes seem to reduce the possibility of discovering Neolithic sites, and thereby also contribute to the uneven site distribution. This is especially true for forest and grassland (Table 2).

Neolithic settlement dynamics
4.3.1 Site frequency. On average, there is a frequency of approx. 43 sites per century for the entire Neolithic. However, strong fluctuations can be observed for the individual periods. The Early Neolithic is characterized by a frequency of ca. 33 sites per 100 years. The transition to the Middle Neolithic is marked by a decline to 5, followed by a massive increase to ca. 60 with the onset of the Late Neolithic. 4.3.2 Site distribution and density. The 1.5 km isolines cover about 90% of all Neolithic sites with exact coordinates. Site concentrations are recognisable in the southernmost part of the Northern German Plain south of Leipzig, along the Weiße Elster river between Zeitz and Gera, and north-east of Gera. Furthermore, the 1.5 km isolines describe areas with a locally increased density of sites in the Central Uplands. These are located in the valley of the Weiße Elster near Plauen, north-east of that town, as well as in the area of the watershed to the river Saale. In addition, there is a loose scattering of sites in the Central Uplands which are located outside the 1.5 or 2.5 km isolines (Fig 8A).
Based on the 1.5 km isolines, the distribution of about 70% of the Early Neolithic sites can be described. An elongated concentration south of Leipzig along the Weiße Elster river valley is recognisable. Further site concentrations can be described with the 1.5 km isolines near Zeitz, Eisenberg and north-east of Gera. These local concentrations are connected by the 2.5 km isoline which covers 87% of the Early Neolithic sites. Far outside these isolines there are isolated sites in the Central Uplands (Fig 8B). The majority of the Middle Neolithic sites is loosely scattered, so that the 1.5 km isoline can only cover about 57% of the sites and the 2.5 km isoline about 77%. A rather elongated concentration along the Weiße Elster valley and a small concentration near Zeitz are described. In addition, there are several isolated sites in a greater distance from the Weiße Elster valley. South of Zeitz there are only very few sites dating to the Middle Neolithic ( Fig 8C). In comparison, the distribution of Late Neolithic sites covers a larger area than the Middle Neolithic. The 1.5 km isolines cover 67% of the Late Neolithic sites and describe several areas with a higher density: a large area along the Weiße Elster valley  [83]. The spatial extent of the (refilled) open-cast mines has been modified according to older topographic maps (S2 Table;    and west of it, smaller areas near Zeitz, Eisenberg, Gera, and one on the watershed between the rivers Weiße Elster and Pleiße. Overall, the distribution of Late Neolithic sites is comparable with the Early Neolithic (Fig 8D).

Site Exploitation Territories (SET) and settlement dynamics.
For all sites, a differentiation of the Neolithic periods is possible based on the smallest SETs observed. For the Early and Late Neolithic sites very small minimal SETs with sizes of about 0.7 km 2 were identified. In contrast, the smallest SETs for Middle Neolithic sites cover areas of 1 km 2 . The same pattern is true for the settlements of the three periods ( Fig 9A).
While the immediate surroundings of Early Neolithic sites are on average located at an elevation of 182 m a.s.l., the SETs of Middle Neolithic sites are characterised by lower elevations (141 m a.s.l.). Late Neolithic SETs are characterised by elevations similar to those determined for the Early Neolithic sites (176 m a.s.l.). A separate evaluation of settlement sites results in comparable values for the Early and Middle Neolithic, i.e. 175 and 143 m a.s.l. (Fig 9B). In contrast, the average elevation within the SETs of Late Neolithic settlements is similar to that of Middle Neolithic settlements (157 m a.s.l.).
With regard to slope, differences between the Early and Late Neolithic on the one hand (ca. 1.6˚), and the Middle Neolithic on the other hand (ca. 0.8˚) can be observed. In addition, the average value for Late Neolithic settlement SETs (1.2˚) is lower than the average for all sites from this period (Fig 9C).
By calculating the average valley depth within the SETs for all sites, the Early and Late Neolithic (ca. 38 and 35 m) can be distinguished from the Middle Neolithic (27 m). However, the average valley depth in the vicinity of Early and Late Neolithic settlements is slightly below the  Table). See Fig 2 for  average for all sites for the respective periods. The opposite is true for the Middle Neolithic (Fig 9D).
For all sites, an increased river distance of about 450 m can be observed with the transition from the Early Neolithic (about 520 m) to the Middle Neolithic (980 m). For the Late Neolithic, the average river distance decreased to about 820 m. The values for Early and Late Neolithic settlements did not differ from those observed for all sites from the respective periods. In contrast, the average river distance within SETs of Middle Neolithic settlements is approx. 100 m lower compared to all sites from that period (Fig 9E).  Table. https://doi.org/10.1371/journal.pone.0265835.g009 The descriptive statistics of the proportion of loess-soils in the SETs did not enable any differentiation between the Neolithic periods ( Fig 9F).

Site distribution with respect to modern land use and archaeological surveys.
Modern land use classes are not tied to archaeological criteria. Therefore, these represent an independent sample with regard to the distribution of archaeological sites in both spatial and chronological terms [100]. The analysis of the CLC dataset showed that urban areas, arable land and open-cast mines are positive filters, since the surface is frequently opened up in these areas by various activities such as construction measures, ploughing or mining. This often leads to the discovery of buried sites ( Table 2). In contrast, forest and grassland were identified as negative filters, since they reduce the likelihood of discovering new sites due to dense vegetation and little to no surface disturbances. These negative filters only dominate in the higher altitudes of the southern study area. It is important to acknowledge the widespread distribution of positive filters especially in the Central Uplands, where only few Neolithic sites have been discovered (Fig 7).
Archaeological surveys led to an increase in Neolithic site densities between Zeitz and Gera (Fig 8). However, this only applies to the Early and Late Neolithic. In contrast, occasional surveys in the Central Uplands did not lead to the discovery of Neolithic sites in those areas. The use of aerial photography did not affect the Neolithic site distributions.

Site distribution with respect to archaeological site visibility and preservation of material remains.
Kappa scores of our RF analyses between 0.53 and 0.59 were found for the depth of Neolithic sites, as well as for the presence/absence and depth of Neolithic pottery with respect to selected topographic, climatic, soil physical and geochemical terrain covariates (S6-S8 Tables). These Kappa scores indicate a substantial agreement between the models and the archaeological classifications [125,126]. This means that the visibility of Neolithic sites and the distribution of material groups as shown in Fig 6 were influenced by geomorphological terrain dynamics and weathering. In this context, the most important terrain covariates were elevation-and subsequently correlated climate parameters such as temperature, frost, precipitation etc.-, slope and the proportion of coarse fragments in soils that largely control/are partly controlled by soil erosion and accumulation processes, and copper distribution that is an unspecific proxy for climatic, geological and pedological factors as well as agricultural activities and land cover [168,169]. The only exception is the distribution of stone tools with regard to their depth. Here, the low Kappa score of 0.21 indicates poor agreement (S9 Table), meaning that the applied terrain covariates hardly describe the distribution of (non-)buried stone tools. Consequently, the preservation conditions for stone tools are similar in the Central Uplands and in the southernmost part of the Northern German Plain, since these were not significantly influenced by the geomorphological terrain dynamics. This finding is important, since stone tools occur at basically every Neolithic site. Hence, the distribution of stone tools can be regarded as a reliable proxy for areas that were inhabited during the Neolithic. Furthermore, our RF results show that despite the general plausibility of the chosen accuracy measurements that were based on general assumptions, the parallel application of several quality measurements is necessary. This is caused by their different specifics according to the data conditions. Here, this is demonstrated by the analysis of the distribution of stone tools with regard to their depth, since error rate, Precision, Recall and the F 1 -measure were in line with all previous measurements and showed a positive agreement, whereas this was opposed by the low Kappa value. The influence of the unbalanced class distributions is evident here in the quality measures (Table 1).
Our finding that the lack of pottery in the Central Uplands (Fig 6B) possibly results from poor local preservation conditions confirms previous speculations for the study area [170][171][172]. However, also subsistence strategies must be taken into account for this phenomenon, as well as for the generally low site densities in the Central Uplands that are a well-known phenomenon in Europe [173]. In contrast to the lowlands, the economic use of low mountain ranges was probably limited to transhumance, which tends to leave few archaeological traces due to its seasonal character and smaller groups of people compared to sedentary farming communities. This mobile way of life meant that only the most necessary equipment could be taken along. Consequently, the pottery carried along also had to be reduced to a minimum [174,175]. Therefore, the absence of Neolithic pottery in the Central Uplands has to be considered with respect to (I) the context of local preservation conditions, and (II) an economic strategy that leaves very little pottery behind and may have been practised in these landscapes. However, further research is be needed to decide which of these two factors has a greater impact on the Neolithic record in the study area.

Neolithic floodplain sites.
Neolithic finds are also known from the sediments of the Weiße Elster floodplain [176]. Originally, these Neolithic sites were most likely located on terrace islands or terrace peninsulas that later became part of the floodplain as we see it today [177,178]. However archaeological, pedological and chronological investigations [179,180] demonstrated that the Neolithic floodplain sites are covered by Subboreal to Sub-Atlantic siltclay overbank fines (floodloam). As a consequence, the site density in these areas may be somewhat higher than current archaeological maps indicate [56]. Similarly, missing knowledge about prehistoric settlement patterns due to thick coverages of prehistoric settlements with younger floodloam are also known from other floodplains [181][182][183].
In a nutshell, based on an assessment of former field surveys, CLC and RF analyses, we can confirm that apart from the distribution of Neolithic pottery and sites in floodplain areas the site distribution in the study area reflects the Neolithic settlement dynamics quite well.

Chronological distribution of Neolithic sites.
The uneven distribution of settlements, burial sites, and stray finds over the different Neolithic periods indicates that their archaeological identification varies between these periods (S3 and S5 Tables). Consequently, these varying likelihoods of identifying different site categories may be one of the driving factors behind fluctuations in site frequencies (Table 3). For example, Middle Neolithic sites might be under-represented because they were identified based on pottery that is rarely decorated. To do so, the shape of largely complete vessels has to be considered for their archaeological dating [184][185][186][187][188][189]. Hence, Middle Neolithic pottery can easily be overlooked as soon as vessel shapes can no longer be derived from pottery fragments [41,56,57]. In contrast, Early and Late Neolithic pottery is decorated with characteristic patterns, and therefore easy to identify [187][188][189]. Furthermore, these latter periods are associated with distinctive stone tools: While shoe-last celts (German: Schuhleistenkeil) are typical for the Early Neolithic, faceted axes (German: Facettenaxt) belong to the most characteristic Late Neolithic artifacts [188,189]. Therefore, the proportion of stray finds is remarkably high for these periods and small for the Middle Neolithic (S3 Table). Moreover, the ratio between settlements and burial sites varies for each period, e.g., more settlements than burial sites are known from the Early and Middle Neolithic while the opposite is true for the Late Neolithic (S3 Table). This is due to the fact that shallow graves were common during the former periods, while the deceased were buried in mounds during the Late Neolithic [89,90,190]. In contrast to burial mounds, shallow graves are hardly visible on the surface and therefore discovered usually by chance [108,135,191,192]. Furthermore, the archaeological identification of Late Neolithic houses is difficult due to their construction techniques [90,99,193]. Results from excavations in Profen [193], Seifartsdorf [194] and Lucka [195,196] indicate that Late Neolithic settlements and burial mounds usually exist in close proximity to each other. Accordingly, Late Neolithic settlements are probably under-represented in the study area, especially in the south where mainly stray finds and burials have been discovered so far. After all, various reasons could have caused the comparatively small number of Middle Neolithic settlements: In addition to a real reduction of settlement activities during that period, this also includes the difficulties in identifying pottery from this period and also the fact that Middle Neolithic settlements were somewhat smaller than settlements from the Early or Late Neolithic. Moreover, the associated features (houses, pits, wells, etc.) are often only loosely scattered. Therefore, Middle Neolithic settlements are most likely to be discovered during excavations of large contiguous areas. In our study area, this was illustrated during excavations at Droßdorf [158,197], Großdalzig [198] and Zauschwitz [199]. Because all the discussed aspects are inherent to the Neolithic periods themselves, these phenomena were observed archaeologically everywhere in Central Germany [33,41,56,57,101,171,[200][201][202][203][204]. Finally, to assess whether the observed reduction of Middle Neolithic compared with Early and Late Neolithic finds was real or not, independent archives of former settlement activity would be necessary for comparisons with our settlement record. However, so far only few pollen records with very limited stratigraphical and chronological resolution exist in this region for the Neolithic period [101,172,[205][206][207][208], so that such comparisons cannot be made so far.

Local settlement dynamics and their regional context
The Neolithic settlement dynamics in the catchment of the river Weiße Elster largely match findings from neighbouring landscapes, including the distinctive decline in Middle Neolithic site frequencies (Table 3; Fig 10). Furthermore, it was repeatedly observed across Central Germany that more elevated zones, which were less favourable for agriculture use were inhabited The Vogtland has not been included, because the respective study does not provide concrete data with respect to the number of Neolithic sites discovered in the area [171]. https://doi.org/10.1371/journal.pone.0265835.t003 during the Early and Late Neolithic, but were largely abandoned during the Middle Neolithic (see references in Fig 10). These settlement dynamics are reflected in the results of our SET analyses. With respect to elevation show that the average elevation values within the SETs dropped at the transition from the Early to the Middle Neolithic. Although no other landscape was analysed using SET and former statistical approaches were based on point coordinates and smaller datasets (see references in Fig 10), our results agree with observations from other studies [33,57,203]. Furthermore, our SET analyses point to greater average distances of Middle and Late Neolithic settlements to rivers compared to the Early Neolithic, which is in line with former observations [33,57,209]. However, site distributions with regard to valley depth were not carried out during former studies so far. Accordingly, there is no information whether the dynamics observed in our study area can also be seen elsewhere. In addition, our SET analyses illustrate that Neolithic settlement dynamics were not linked with different proportions of soils on loess. This is in line with observations from adjacent study areas in Central Germany [57,101].
Another common feature with other regional studies is the concentration of the highest site densities mostly along the floodplains of large rivers (Fig 8; S1 Fig) [23, 210]. This phenomenon has been discussed against the background of the distribution of forests and open land, although there are no pollen records available that can be used to reconstruct the Neolithic vegetation coverage in Central Germany in detail [56]. Hypothetical mappings of forest and open land have been prepared using medieval records, onomastic data and soil maps. Based on these studies, a forest coverage is assumed for the Neolithic floodplains [210][211][212][213]. In our study area, this is supported by tree logs especially of Quercus, but also Fraxinus, Alnus, Pinus and Ulmus with Neolithic ages in floodplain sediments of the lower Weiße Elster valley [214,215]. In the northern part of our study area the forest is supposed to haven been contrasted with the adjacent loess areas, for which a thin and patchy forest vegetation is assumed [56, 97,210,216]. Therefore, the margins between the floodplains and the neighbouring loess-covered slopes may have been preferred settlement areas due to economic benefits, since they offered easy access to both fertile soils and adjacent forests. The latter could be exploited for the extraction of wood for the construction of houses and wells, or be used as fuel, hunting grounds or for forest grazing by cattle, sheep, goats and pigs. The rivers themselves might have been used for fishing and fresh water supply [23,210,217,218]. Especially in the northern part of the study area fresh water resources were probably rare outside the main river valleys. This could have been caused by the relatively flat topography, leading to an extensive lack of small tributary valleys that potentially offer surface water resources by small creeks and sources. Furthermore, the climate of central Europe before ca. 5 ka was generally drier compared with today [5], so that several small current rivers were dry during that period [106]. Altogether, this suggests that these multi-purpose near-river locations were probably preferred settlement areas to compensate economic stresses in these environments (Fig 8) [210]. Recently, this traditional explanation was somehow challenged by large-scale excavations in our study area, during which burial sites and settlements were discovered in the watershed areas between the rivers Weiße Elster and Saale [193] and between the rivers Weiße Elster and Pleiße [219,220]. However, Neolithic wells were discovered in both, the large river valleys and those regions far away from the large valleys. In the latter areas, these, could have helped to provide water resources for the Neolithic settlers and therefore facilitated settlement activities [220,221].

Driving factors of the local settlement dynamics
The main drivers of Neolithic settlement dynamics are still being debated. Among other things, it has been hypothesised that strong population growth during the Early and Late Neolithic may have contributed to the settlement on lower-yielding soils in the Central Uplands [170,201,202,222]. In addition, a connection with climatic conditions has been considered: A generally drier climate during the Early and Late Neolithic could have been compensated by settlement activities in the more humid Central Uplands and well constructions in the southernmost part of the Northern German Plain [223,224]. In fact, most Neolithic wells in the study area have been dated to the Early or Late Neolithic using dendrochronology and radiocarbon dating [158]. Finally, the concentration of Middle Neolithic settlement on the lowlands in the southernmost part of the Northern German Plain was interpreted as an indication that the comparatively high annual temperatures and corresponding longer vegetation periods in this landscape were of special importance for the applied subsistence strategies [99].

Conclusion
Based on a dataset of 1365 Neolithic sites that was compiled in the frame of a joint geoarchaeological research project in Central Germany, we combined methods of archaeological source criticism with machine learning to discuss the large-scale Neolithic settlement dynamics in the immediate catchment of the Weiße Elster river in Central Germany. Our research lead to the following conclusions: • The use of local area files (German: Ortsakten) proved to be extremely valuable. It turned out that only half of the actually existing sites have been published so far.
• Analyses of site distributions with respect to modern land use and field surveys indicate that the basic trends of Neolithic settlement dynamics can be derived from the recorded sites. However, by investigating the influence of terrain covariates on the depth of sites using the Random Forests algorithm, we demonstrate that the visibility of Neolithic sites is influenced by geomorphological terrain dynamics: While the majority of the sites is buried in the southernmost part of the Northern German Plain, the opposite is true for the Central Uplands.
Furthermore, Random Forests algorithm is able to reproduce the distribution of sites with (-out) pottery. These results can be explained in different ways: Firstly, the preservation of pottery in the Central Uplands may be not as good as in the southernmost part of the Northern German Plain. Secondly, the low amount of Neolithic pottery in the Central Uplands may also be a result of economic strategies that leave very little pottery behind in general. Nevertheless, the Random Forests algorithm offers new opportunities for archaeological source criticism, which is crucial for identifying potential biases in large archaeological datasets and deriving settlement dynamics.
• The distribution of Early Neolithic sites in the lowlands of the southernmost part of the Northern German Plain and the mountainous landscape between Zeitz and Gera is characterized by clusters along the Weiße Elster River. In contrast, the geographic focus of Middle Neolithic land use was restricted to the northern lowlands, whereas the Late Neolithic site distribution pattern resembled that of the Early Neolithic. The site frequencies are marked by a decline in the Middle Neolithic, and a re-increase during the Late Neolithic. These observed trends could at least partly be consequences of diachronically different properties of the archaeological record.
• In the southernmost part of the Northern German Plain, Neolithic land use generally took place at the transition between the floodplain of the Weiße Elster River and loess-covered soils outside the valley. These settlement patterns are common in Central Germany, and were explained with environmental economic factors. However, both settlements and burial sites have been discovered recently in the watershed areas between Weiße Elster and Saale River as well as between the rivers Weiße Elster and Pleiße. Accordingly, the low site densities in greater distances from rivers might be a reflection of a lack of research.
• Analyses of SETs essentially reflect the settlement dynamics between the southernmost part of the Northern German Plain and the Central Uplands. Early Neolithic settlements were located closer to river valleys, and were situated in relatively elevated areas above the river valleys. In contrast, settlements dating to later periods were not located as high above the river valleys, but in larger distances from the latter. However, Neolithic settlement dynamics were not linked with different proportions of soils on loess.
• We identified potential for future research in the region: Numerous Neolithic sites in the surroundings of Gera could not be assigned to any period so far, but might be re-dated by means of archaeological methods. Moreover, the Late Neolithic in the Central Uplands is mainly known through burials and stray finds.
Generally, our study demonstrates the high value of systematic studies of diachronic archaeological settlement patterns to understand varying regional human activity in space and time. However, since such records can show different properties and sensitivities through time, complementary studies of other environmental or historic archives documenting regional human impact could help to identify possible intrinsic biases. In the future complementary environmental archives will be analysed to get deeper insights into early settlement processes and their environmental background. A key question than will be to which degree human activities and climatic factors were responsible for fluvial soil deposits.
Supporting information S1  Table. Proportion of (non-)intentional as well as unknown modes of discovery associated with each Neolithic period. Table. Results of random forests analyses with respect to interrelation of terrain covariates and site depth. Table. Results of random forests analyses with respect to interrelation of terrain covariates and the distribution of material groups. (XLSX) S8 Table. Results of random forests analyses with respect to interrelation of terrain covariates and the distribution of (non-)buried pottery. Table. Results of random forests analyses with respect to interrelation of terrain covariates and the distribution of (non-)buried stone tools.